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Abstract. We present a simple theory for the evolution of initially compact clusters in a tidal 
field. The fundamental ingredient of the model is that a cluster conducts a constant fraction 
of its own energy through the half-mass radius by two-body interactions every half-mass 
relaxation time. This energy is produced in a self-regulative way in the core by an (unspeci- 
fied) energy source. We find that the half-mass radius increases during the first part (roughly 
. half) of the evolution and decreases in the second half, while the escape rate is constant and 

' set by the tidal field. We present evolutionary tracks and isochrones for clusters in terms 

\ of cluster half-mass density, cluster mass and galacto-centric radius. We find substantial 

agreement between model isochrones and Milky Way globular cluster parameters, which 
suggests that there is a balance between the flow of energy and the central energy produc- 
tion for almost all globular clusters. We also find that the majority of the globular clusters 
are still expanding towards their tidal radius. Finally, a fast code for cluster evolution is 
?—j . presented. 
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£> . 1. Introduction tracks and isochrones for clusters (i.e. not for 

!■/"") ' the stars within them!) evolving in a tidal field. 

Capturing cluster evolution in equations is This foms the theoretical f rame work to ex- 
complex because several processes, including plain empirically established correlations be- 
lO ■ two-body relaxation, interactions with binary tween stmctural p arame ters and their environ- 
stars, escape across the tidal boundary, and the ment afj found for MOky Way globular dusters 
internal evolution and mass-loss of single and (e g Djorgovski 19 9 5; McLaughlin 2000) and 
m ; binary stars are all at work at the same time. for extl -a-galactic globular cluster systems (e.g. 

However, it is desirable to have a simple pa- Jordan et al 2005; McLaughlin et al . 2 008; 

rameterisation of the evolution of some funda- Harris et al 2010) We do not aim tQ explain 

mental cluster parameters such as mass and ra- the shape and dependence on environment of 

dms (e.g. Pneto & Gnedin 2008). the globular cluster mass function . 

Here we provide a physically motivated 
and simple prescription for the behaviour of The other principal exclusion is the evolu- 

fhe half-mass radius and tidal radius (i.e. tion of the core parameters (i.e. the core mass 

mass). With this we construct evolutionary and radius). Although the core is the place 

where the energy is produced, the assumption 

Send offprint requests to: M. Gieles@surrey.ac.uk of excluding it in the model is justified by 
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the discovery of Henon (1975) that the rate of 
flow of energy is controlled by the system as 
a whole, and not by the core. In Henon's pic- 
ture the mechanism of energy generation in the 
core is self-regulatory and so we can assume 
that the core produces the right amount of en- 
ergy required by the system as a whole. This is 
comparable to the self-regulative energy pro- 
duction in stellar cores, which was first realised 
by Eddington. The application of this idea to 
stellar dynamics was a breakthrough allowing 
modellers to overcome the core collapse phase. 

From A^-body simulations of the long term 
(post-collapse) evolution of single-mass clus- 
ters, i.e. where stars have the same mass, it was 
found that binary stars act as the energy source 
(Giersz & Heggie 1994; Baumgardt et al. 
2002). In models with more realistic initial 
conditions these binaries are usually consid- 
ered to be primordial. Other mechanisms of en- 
ergy generation have been considered, includ- 
ing the action of a central intermediate-mass 
black hole (Baumgardt et al. 2004; Heggie 
et al. 2007) 

Mass-loss from stellar evolution can also 
provide the energy for the dynamical evolu- 
tion of clusters (Gieles et al. 2010). Typically 
most of the mass is lost from the most mas- 
sive stars in the cluster, that reside in the cluster 
core as the result of mass segregation. The re- 
sulting energy production works together with 
binaries in driving an expansion on relaxation 
time-scales (Gieles et al. 2010). Indeed, it may 
even dominate, as in the specific example of 47 
Tuc, a high-concentration cluster in which the 
evolution of the core- and half-mass radii ap- 
pear to be little affected by the primordial bi- 
nary population (Giersz & Heggie 201 1). 

These self-regulatory mechanisms of en- 
ergy generation take time to establish the bal- 
ance between the energy generated in the core 
and the energy requirements of the overall evo- 
lution of the cluster. We refer to the subsequent 
evolution as 'balanced'. In much research on 
cluster dynamics this kind of evolution is usu- 
ally associated with 'post-collapse' evolution, 
but this term is ambiguous; the phrase 'post- 
collapse' might be used for the entire evolution 
after the end of mass segregation of massive 
stars (> 10 M ), whilst others use it to refer 



to the evolution following the decrease in the 
core radius after several Gyrs. For this reason 
we prefer the term 'balanced evolution' . 

In Section 2 we present the model that uni- 
fies two models of Henon. The evolution in 
the initial phase resembles that of the model 
for the isolated cluster which expands with lit- 
tle loss of stars (Henon 1961, hereafter H61) 
and near final dissolution the cluster resembles 
the 'homologous cluster' (Henon 1965, here- 
after H65), which contracts at a constant den- 
sity. Here the model is described in its sim- 
plest form and we refer to (Gieles, Heggie, & 
Zhao 2011, hereafter Gil) for more details. 
In Section 3 the model is compared to param- 
eters of the Milky Way globular clusters. In 
Section 4 a fast code for cluster evolution based 
on the theory presented in this work is pre- 
sented. 

2. Model 

Our attempt to unify the evolution of mass and 
radius of the two models of Henon begins with 
one of the physical properties which the two 
models have in common, i.e. a flux of energy 
at the half-mass radius which is fed by an en- 
ergy source in the core. We constrain ourselves 
to the energy flow at this radius, because we 
can then construct a relatively simple set of re- 
lations for the behaviour of the bulk properties 
of the cluster. We adopt Henon's idealisation of 
systems in which all stars have the same mass 
m. We estimate the total energy of the cluster 
as usual by 



where M is the mass of the cluster, rj, is the 
half-mass radius, and a ^ 0.2 is a 'form fac- 
tor'. 

In constructing a unified approximate 
model which includes the transition from 
nearly isolated evolution to tidally limited evo- 
lution, we assume that in both phases there is 
an energy flow due to two-body relaxation of 
magnitude 

E t 

Tm = ■ (2) 
\E\ T rh 
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Here r r h is the half-mass relaxation time-scale 
and £ is a constant that can be interpreted as 
the efficiency of energy conduction. In Henon's 
models £ 0.08 1 and from numerical simula- 
tions Alexander & Gieles (2012) find f 0.1. 
We approximate the expression for t± by as- 
suming that the Coulomb logarithm is con- 
stant, such that (Spitzer 1987) 

T rh oc Nt ci . (3) 

Here N = M/m is the total number of stars and 
r cr is the crossing time of stars in the cluster at 
the half-mass radius. We define r cr as 

Tcr = (Gph)- 1/2 , (4) 

with ph = 3M/ (inrfy the cluster density within 
the half-mass radius. 

Before proceeding further, we shall change 
the variables in which the total energy E is ex- 
pressed, because this will facilitate the further 
development of our model. Instead of using the 
half-mass radius, r^, we shall use T cr such that 
E oc -M 5/3 t ct 2/3 . Putting this together with our 
assumption about the energy flux (equation 2) 
we find 

_5M + 2t C£ = jL 
3 M 3 r cr r rh 

We recall that r r h oc Mr cr and thus equa- 
tion (5) has two variables: M and r cr . The dif- 
ferential equation can be solved by relating M 
to r cr . Because M depends only on the orbit 
and is, to good approximation, independent of 
the cluster mass and radius (Lee & Ostriker 
1987; Gieles & Baumgardt 2008) we can write 
for the dimensionless escape rate (M/M)T T h = 
-(3/5)^T cr /T crl . Here T cr i is the maximum r cr 
which depends on the tidal field: if the tides 
are weak, the cluster can expand to larger r cr . 
Combining this dimensionless escape rate with 
equation (5) we find the dimensionless expan- 
sion rate (t cr /T cr )r rh = (3/2)f(l - T cr /T crl ). 
Dividing M by t cr we find the surprisingly sim- 
ple relation between M and T cr 

AM 2 M 

j — = 7 ■ (°) 

ar cr j r cr - T cr i 



Integration gives an expression for r cr (M, r cr i), 
i.e. the isochrones, and this, combined with 
equation (5), can be used to get the time- 
dependent solutions M{f) and r cr (f), i.e. the 
evolutionary tracks. We do not give the func- 
tional forms here, but instead refer the reader 
to Gl 1. In the next section we proceed with a 
direct comparison between cluster isochrones 
and Milky Way globular clusters parameters. 

3. Milky Way globular clusters 

3.1. Are clusters still expanding? 

In order to see if these types of predictions are 
relevant for real globular clusters we compare 
our results to the globular clusters of the Milky 
Way. We use the 2003 version of the Harris 
(1996) catalogue which contains entries for 
150 globular clusters, and for 141 of them a lu- 
minosity, half-light radius and galacto-centric 
radius determination are available. To convert 
luminosity to mass we adopt a mass-to-light 
ratio of 2 (McLaughlin & van der Marel 2005) 
and we multiply the projected half-light radius 
by 4/3 to correct for the effect of projection 
(Spitzer 1987) and get an estimate for r^. 

The first thing we determine from the data 
is the fraction of globular clusters that are in the 
expansion dominated phase. We define the end 
of the expansion phase as the moment where 
the time derivative of T rri is zero, which is 
when the fractional change in M has the same 
magnitude as the fractional change in r cr , i.e. 
dlnM/dlnr cr = -1 (equation 3). Combined 
with equation (6) we then find that this hap- 
pens when M/M = (2/7) 2/5 0.6. In time 
this is when the cluster has evolved for 40% of 
its total life, because of the linear decrease of 
M. The end of the expansion phase therefore 
depends on M and the mass-loss rate M, which 
depends on the orbit. 

We express M in terms of the galactocen- 
tric radius Rq. We assume that the Milky Way 
halo is an isothermal sphere and approximate 
the Rq dependent mass-loss rate by 

MR C -20 M Myr" 1 kpc, (7) 



1 In Gl 1 the values for the isolated and homolo- 
gous cluster are derived from Henon's papers. 



which is in reasonable agreement with the 
evaporation rates found in both N-body and 
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Fig. 1. Isochrones based on the balanced evolution model of Section 2. In the left panel 5 
isochrones for different Rq values are shown and in the right panel 5 isochrones for different 
M. Both panels show the 141 globular clusters in the Harris (1996) catalogue from which M, rj, 
and Rq are available. In both panels different symbols are used for data in different Rq and M 
regimes. 



Fokker-Planck models with a globular cluster 
type stellar mass function (Gil). We adopt an 
age of 13 Gyr for all clusters. Clusters that are 
now at the end of the expansion phase have 
60% of their evolution, or 1.5 x 13 Gyr, ahead 
of them. The remaining life-time, or life ex- 
pectancy, is defined as M/M such that with 
equation (7) we find that clusters with a mass 

M>lO 5 M (!|p) (8) 

are still in the expansion-dominated phase. 
This relation is satisfied by 93 of the 141 clus- 
ters (i.e. roughly 2/3). It follows that the re- 
maining 48 clusters (roughly 1/3) have ex- 
panded to the tidal boundary and are in the 
evaporation-dominated phase. This perhaps 
surprising result has some interesting conse- 
quences. The most important one is that the 
present day densities of the majority of the 
globular clusters follow (roughly) from the 
self-similar expansion model for isolated clus- 
ters: r r h oc Mt c1 - si constant, i.e. ph 00 M 2 . 
This scaling relation is caused by internal two- 



body relaxation and is independent of the tidal 
field; therefore a similar scaling, with the same 
proportionality, should also hold for extra- 
galactic clusters. Moreover, in extra-galactic 
cluster samples the fraction of clusters in the 
expansion-dominated phase is probably larger; 
they are easier to detect because they have (on 
average) higher mass (equation 8). 

1 /2 

The prediction that a p h ' oc M scaling 
must hold for the majority of the Milky Way 
globular clusters is one of the main results of 
this work. 

3.2. Isochrones 

Because all globular clusters have roughly the 
same age we focus on isochrones with an age 
of 13 Gyr, rather than the evolutionary tracks. 
For a given age the isochrones can be ex- 
pressed as ph(M, Rq) and the detailed results 
are given in the appendices of Gil. In Fig. 1 
we show the isochrones ph(M) for several val- 
ues of Rq (left) and ph(^G) for several values of 
M (right) diagrams together with the 141 glob- 
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ular clusters for which data are available in the 
Harris catalogue. 

In the left panel isochrones for clusters at 
different Rq between lkpc and lOOkpc are 
shown. The isochrones roughly encompass the 
data. The lOOkpc isochrone clearly shows the 

1 /2 

asymptotic p h oc M behaviour following from 
expansion, which roughly follows the lower 
envelope of data points. In the outer halo the 
tidal field is so weak that all clusters with 
M > 10 4 M o are still expanding towards their 
tidal boundary. In the right panel the densities 
are shown as a function of R G together with 
five isochrones for different masses. These 
isochrones also roughly encompass the data. 
The asymptotic behaviour of the isochrones in 
both diagrams is given by labels in the two di- 
agrams. 

4. A fast code for cluster evolution 

The model presented here makes several ap- 
proximations in order to facilitate simple an- 
alytical results. An improved version of the 
model, which includes the escape of stars in 
the isolated regime, the small N dependence 
in the Coulomb logarithm and the delayed es- 
cape of stars due to the anisotropic tidal field 
(Fukushige & Heggie 2000) is presented in 
Alexander & Gieles (2012). Here the differ- 
ential equations for M and n, are solved nu- 
merically with a Runge-Kutta solver, which 
gives near instantaneous results for M(t) and 
rh(f) as a function of the tidal field strength 
and the initial cluster parameters which can be 
specified on the command line. The code ac- 
curately reproduces the results of iV-body in- 
tegrations of single-mass clusters (Alexander 
& Gieles 2012) and is publicly available on 
https://github.com/emacss/emacss. Future ver- 
sions of the code will be able to reproduce the 
evolution of more realistic stellar clusters, in- 
cluding a stellar mass function and the effects 
of stellar evolution. 
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